On this page

Spectral Clustering

Spectral Clustering

Overview

The current notebook will demonstrate a simplified machine learning approach to observe the change in a lake water’s extent across time. In order to identify the water, we can use spectral clustering to classify each grid cell into a category based on the similarity of the combined set of pixels across wavelength-bands in our image stacks.

Our example approach uses a version of spectral clustering from dask_ml that is a scalable equivalent of what is available in scikit-learn. We will begin this approach with a single image stack and then conduct a direct comparison on the results from different time points.

This workflow uses data from Microsoft Planetary Computer but it can be adapted to work with any data ingestion approach from this cookbook.

Prerequisites

Concepts

Importance

Notes

Data Ingestion - Planetary Computer

Necessary

scikit-learn

Helpful

Spectral clustering

dask_ml

Helpful

Spectral clustering at scale

  • Time to learn: 20 minutes.

Imports

# Data
import numpy as np
import odc.stac
import pandas as pd
import planetary_computer
import pystac_client
import xarray as xr
from dask.distributed import Client
from pystac.extensions.eo import EOExtension as eo

# Analysis
from dask_ml.cluster import SpectralClustering

# Viz
import hvplot.xarray

/home/runner/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/dask/dataframe/_pyarrow_compat.py:17: FutureWarning: Minimal version of pyarrow will soon be increased to 14.0.1. You are using 12.0.1. Please consider upgrading.
  warnings.warn(

Loading Data

Let’s start by loading some Landsat data. These steps are covered in the Data Ingestion - Planetary Computer prerequisite.

Search the catalog

catalog = pystac_client.Client.open(
    "https://planetarycomputer.microsoft.com/api/stac/v1",
    modifier=planetary_computer.sign_inplace,
)

bbox = [-118.89, 38.54, -118.57, 38.84]  # Region over a lake in Nevada, USA
datetime = "2017-06-01/2017-09-30"  # Summer months of 2017
collection = "landsat-c2-l2"
platform = "landsat-8"
cloudy_less_than = 1  # percent

search = catalog.search(
    collections=["landsat-c2-l2"],
    bbox=bbox,
    datetime=datetime,
    query={"eo:cloud_cover": {"lt": cloudy_less_than}, "platform": {"in": [platform]}},
)
items = search.get_all_items()
print(f"Returned {len(items)} Items:")
[[i, item.id] for i, item in enumerate(items)]

/home/runner/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/pystac_client/item_search.py:849: FutureWarning: get_all_items() is deprecated, use item_collection() instead.
  warnings.warn(

Returned 3 Items:

[[0, 'LC08_L2SP_042033_20170718_02_T1'],
 [1, 'LC08_L2SP_042033_20170702_02_T1'],
 [2, 'LC08_L2SP_042033_20170616_02_T1']]

Load a dataset

item = items[1]  # select one of the results

assets = []
for _, asset in item.assets.items():
    try:
        assets.append(asset.extra_fields["eo:bands"][0])
    except:
        pass

cols_ordered = [
    "common_name",
    "description",
    "name",
    "center_wavelength",
    "full_width_half_max",
]
bands = pd.DataFrame.from_dict(assets)[cols_ordered]
bands
common_name description name center_wavelength full_width_half_max
0 red Visible red OLI_B4 0.65 0.04
1 blue Visible blue OLI_B2 0.48 0.06
2 green Visible green OLI_B3 0.56 0.06
3 nir08 Near infrared OLI_B5 0.87 0.03
4 lwir11 Long-wave infrared TIRS_B10 10.90 0.59
5 swir16 Short-wave infrared OLI_B6 1.61 0.09
6 swir22 Short-wave infrared OLI_B7 2.20 0.19
7 coastal Coastal/Aerosol OLI_B1 0.44 0.02 
ds_2017 = odc.stac.stac_load(
    [item],
    bands=bands.common_name.values,
    bbox=bbox,
    chunks={},  # <-- use Dask
).isel(time=0)

Retain CRS Attribute

epsg = item.properties["proj:epsg"]
ds_2017.attrs["crs"] = f"epsg:{epsg}"

da_2017 = ds_2017.to_array(dim="band")
da_2017

<xarray.DataArray (band: 8, y: 1128, x: 950)>
dask.array<stack, shape=(8, 1128, 950), dtype=uint16, chunksize=(1, 1128, 950), chunktype=numpy.ndarray>
Coordinates:
  * y            (y) float64 4.301e+06 4.301e+06 ... 4.267e+06 4.267e+06
  * x            (x) float64 3.353e+05 3.353e+05 ... 3.637e+05 3.638e+05
    spatial_ref  int32 32611
    time         datetime64[ns] 2017-07-02T18:33:06.200763
  * band         (band) object 'red' 'blue' 'green' ... 'swir22' 'coastal'
Attributes:
    crs:      epsg:32611

Reshaping Data

The shape of our data is currently n_bands, n_y, n_x. In order for dask-ml / scikit-learn to consume our data, we’ll need to reshape our image stacks into n_samples, n_features, where n_features is the number of wavelength-bands and n_samples is the total number of pixels in each wavelength-band image. Essentially, we’ll be creating a vector of pixels out of each image, where each pixel has multiple features (bands), but the ordering of the pixels is no longer relevant to the computation.

By using xarray methods to flatten the data, we can keep track of the coordinate labels ‘x’ and ‘y’ along the way. This means that we have the ability to reshape back to our original array at any time with no information loss!

flattened_xda = da_2017.stack(z=("x", "y"))  # flatten each band
flattened_t_xda = flattened_xda.transpose("z", "band")
flattened_t_xda

<xarray.DataArray (z: 1071600, band: 8)>
dask.array<transpose, shape=(1071600, 8), dtype=uint16, chunksize=(1071600, 1), chunktype=numpy.ndarray>
Coordinates:
    spatial_ref  int32 32611
    time         datetime64[ns] 2017-07-02T18:33:06.200763
  * band         (band) object 'red' 'blue' 'green' ... 'swir22' 'coastal'
  * z            (z) object MultiIndex
  * x            (z) float64 3.353e+05 3.353e+05 ... 3.638e+05 3.638e+05
  * y            (z) float64 4.301e+06 4.301e+06 ... 4.267e+06 4.267e+06
Attributes:
    crs:      epsg:32611

Standardize Data

Now that we have the data in the correct shape, let’s standardize (or rescale) the values of the data. We do this to get all the flattened image vectors onto a common scale while preserving the differences in the ranges of values. Again, we’ll demonstrate doing this first in NumPy and then xarray.

with xr.set_options(keep_attrs=True):
    rescaled_xda = (flattened_t_xda - flattened_t_xda.mean()) / flattened_t_xda.std()
rescaled_xda

<xarray.DataArray (z: 1071600, band: 8)>
dask.array<truediv, shape=(1071600, 8), dtype=float64, chunksize=(1071600, 1), chunktype=numpy.ndarray>
Coordinates:
    spatial_ref  int32 32611
    time         datetime64[ns] 2017-07-02T18:33:06.200763
  * band         (band) object 'red' 'blue' 'green' ... 'swir22' 'coastal'
  * z            (z) object MultiIndex
  * x            (z) float64 3.353e+05 3.353e+05 ... 3.638e+05 3.638e+05
  * y            (z) float64 4.301e+06 4.301e+06 ... 4.267e+06 4.267e+06
Attributes:
    crs:      epsg:32611

Info

Above, we are using a context manager “with xr.set_options(keep_attrs=True):” to retain the array’s attributes through the operations. That is, we want any metadata like ‘crs’ to stay with our result so we can use ‘geo=True’ in our plotting.

As rescaled_xda is still a Dask object, if we wanted to actually run the rescaling at this point (provided that all the data can fit into memory), we would use rescaled_xda.compute().

ML pipeline

Now that our data is in the proper shape and value range, we are ready to conduct spectral clustering. Here we will use a version of spectral clustering from dask_ml that is a scalable equivalent to operations from Scikit-learn that cluster pixels based on similarity (across all wavelength-bands, which makes it spectral clustering by spectra!)

client = Client(processes=False)
client

Client

Client-e85a32dd-e6e5-11ee-8e38-6045bde8e8e2

Connection method: Cluster object Cluster type: distributed.LocalCluster
Dashboard: http://10.1.1.51:8787/status

Cluster Info

Now we will compute and persist the rescaled data to feed into the ML pipeline. Notice that our X matrix below has the shape: n_samples, n_features as discussed earlier.

X = client.persist(rescaled_xda)
X.shape

(1071600, 8)

First we will set up the model with the number of clusters, and other options.

clf = SpectralClustering(
    n_clusters=4,
    random_state=0,
    gamma=None,
    kmeans_params={"init_max_iter": 5},
    persist_embedding=True,
)

This next step is the slow part. We’ll fit the model to our matrix X. Depending on your setup, it could take seconds to minutes to run depending on the size of our data.

%time clf.fit(X)

/home/runner/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/distributed/client.py:3169: UserWarning: Sending large graph of size 81.80 MiB.
This may cause some slowdown.
Consider scattering data ahead of time and using futures.
  warnings.warn(

/home/runner/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/dask/base.py:1462: UserWarning: Running on a single-machine scheduler when a distributed client is active might lead to unexpected results.
  warnings.warn(

CPU times: user 22.4 s, sys: 13.1 s, total: 35.5 s
Wall time: 29.3 s

SpectralClustering(gamma=None, kmeans_params={'init_max_iter': 5}, n_clusters=4,
                   persist_embedding=True, random_state=0)
In a Jupyter environment, please rerun this cell to show the HTML representation or trust the notebook.
On GitHub, the HTML representation is unable to render, please try loading this page with nbviewer.org.

Let’s check the shape of the result:

labels = clf.assign_labels_.labels_.compute()
labels.shape

(1071600,)

labels

array([0, 2, 2, ..., 2, 2, 2], dtype=int32)

The result is a single vector of cluster labels.

Un-flattening

Once the computation is done, we can use the coordinates of our input array to restack our output array back into an image. Again, one of the main benefits of using xarray for this stacking and unstacking is that it keeps track of the coordinate information for us.

Since the original array is n_samples by n_features (90000, 6) and the cluster label output is (90000,), we just need the coordinates from one of the original features in the shape of n_samples. We can just copy the coordinates from the first input feature and populate is with our output data:

template = flattened_t_xda[:, 0]
output_array = template.copy(data=labels)
output_array

<xarray.DataArray (z: 1071600)>
array([0, 2, 2, ..., 2, 2, 2], dtype=int32)
Coordinates:
    spatial_ref  int32 32611
    time         datetime64[ns] 2017-07-02T18:33:06.200763
    band         <U3 'red'
  * z            (z) object MultiIndex
  * x            (z) float64 3.353e+05 3.353e+05 ... 3.638e+05 3.638e+05
  * y            (z) float64 4.301e+06 4.301e+06 ... 4.267e+06 4.267e+06
Attributes:
    crs:      epsg:32611

With this new output array with coordinates copied from the input array, we can unstack back to the original x and y image dimensions by just using .unstack().

unstacked_2017 = output_array.unstack()
unstacked_2017

/home/runner/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/numpy/core/numeric.py:407: RuntimeWarning: invalid value encountered in cast
  multiarray.copyto(res, fill_value, casting='unsafe')

<xarray.DataArray (x: 950, y: 1128)>
array([[0, 2, 2, ..., 0, 0, 0],
       [0, 0, 0, ..., 0, 0, 0],
       [2, 0, 2, ..., 0, 0, 0],
       ...,
       [2, 2, 2, ..., 2, 2, 2],
       [2, 2, 2, ..., 2, 2, 2],
       [2, 2, 0, ..., 2, 2, 2]], dtype=int32)
Coordinates:
  * x            (x) float64 3.353e+05 3.353e+05 ... 3.637e+05 3.638e+05
  * y            (y) float64 4.301e+06 4.301e+06 ... 4.267e+06 4.267e+06
    spatial_ref  int32 32611
    time         datetime64[ns] 2017-07-02T18:33:06.200763
    band         <U3 'red'
Attributes:
    crs:      epsg:32611

Finally, we can visualize the results! By hovering over the resulting imge, we can see that the lake water has been clustered with a certain label or ‘value’.

raw_plot_2017 = da_2017.sel(band="red").hvplot.image(
    x="x", y="y", geo=True, xlabel="lon", ylabel="lat", datashade=True, cmap="greys", title="Raw Image 2017",
)

result_plot_2017 = unstacked_2017.hvplot(
    x="x", y="y", cmap="Set3", geo=True, xlabel="lon", ylabel="lat", colorbar=False,  title="Spectral Clustering 2017",
)

raw_plot_2017 + result_plot_2017

Spectral Clustering for 1988

We have conducted the spectral clustering for 2017 and now we want to compare this result to the lake in 1988. Let’s load data from 1988 and run the same analysis as above.

We will use the same catalog, but we will search it for a different point in time and different Landsat mission

Load the data

bbox = [-118.89, 38.54, -118.57, 38.84]  # Region over a lake in Nevada, USA
datetime = "1988-06-01/1988-09-30"  # Summer months of 1988
collection = "landsat-c2-l2"
platform = "landsat-5"  # Searching through an earlier landsat mission
cloudy_less_than = 1  # percent

search = catalog.search(
    collections=["landsat-c2-l2"],
    bbox=bbox,
    datetime=datetime,
    query={"eo:cloud_cover": {"lt": cloudy_less_than}, "platform": {"in": [platform]}},
)

items = search.get_all_items()
item = items[1]  # select one of the results

/home/runner/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/pystac_client/item_search.py:849: FutureWarning: get_all_items() is deprecated, use item_collection() instead.
  warnings.warn(

Notice that Landsat 5 data from 1988 has slightly different spectra than Landsat 8 from 2017. Details like this are important to keep in mind when performing analyses that directly compare across missions.

assets = []
for _, asset in item.assets.items():
    try:
        assets.append(asset.extra_fields["eo:bands"][0])
    except:
        pass

cols_ordered = [
    "common_name",
    "description",
    "name",
    "center_wavelength",
    "full_width_half_max",
]
bands = pd.DataFrame.from_dict(assets)[cols_ordered]
bands
common_name description name center_wavelength full_width_half_max
0 red Visible red TM_B3 0.66 0.06
1 blue Visible blue TM_B1 0.49 0.07
2 lwir Long-wave infrared TM_B6 11.45 2.10
3 green Visible green TM_B2 0.56 0.08
4 nir08 Near infrared TM_B4 0.83 0.14
5 swir16 Short-wave infrared TM_B5 1.65 0.20
6 swir22 Short-wave infrared TM_B7 2.22 0.27 
ds_1988 = odc.stac.stac_load(
    [item],
    bands=bands.common_name.values,
    bbox=bbox,
    chunks={},  # <-- use Dask
).isel(time=0)

epsg = item.properties["proj:epsg"]
ds_1988.attrs["crs"] = f"epsg:{epsg}"

da_1988 = ds_1988.to_array(dim="band")
da_1988

<xarray.DataArray (band: 7, y: 1128, x: 950)>
dask.array<stack, shape=(7, 1128, 950), dtype=uint16, chunksize=(1, 1128, 950), chunktype=numpy.ndarray>
Coordinates:
  * y            (y) float64 4.301e+06 4.301e+06 ... 4.267e+06 4.267e+06
  * x            (x) float64 3.353e+05 3.353e+05 ... 3.637e+05 3.638e+05
    spatial_ref  int32 32611
    time         datetime64[ns] 1988-07-02T18:03:59.010013
  * band         (band) object 'red' 'blue' 'lwir' ... 'nir08' 'swir16' 'swir22'
Attributes:
    crs:      epsg:32611

Reshape and Standardize

flattened_xda = da_1988.stack(z=("x", "y"))
flattened_t_xda = flattened_xda.transpose("z", "band")
with xr.set_options(keep_attrs=True):
    rescaled_xda = (flattened_t_xda - flattened_t_xda.mean()) / flattened_t_xda.std()
rescaled_xda

<xarray.DataArray (z: 1071600, band: 7)>
dask.array<truediv, shape=(1071600, 7), dtype=float64, chunksize=(1071600, 1), chunktype=numpy.ndarray>
Coordinates:
    spatial_ref  int32 32611
    time         datetime64[ns] 1988-07-02T18:03:59.010013
  * band         (band) object 'red' 'blue' 'lwir' ... 'nir08' 'swir16' 'swir22'
  * z            (z) object MultiIndex
  * x            (z) float64 3.353e+05 3.353e+05 ... 3.638e+05 3.638e+05
  * y            (z) float64 4.301e+06 4.301e+06 ... 4.267e+06 4.267e+06
Attributes:
    crs:      epsg:32611

Spectral Clustering

X = client.persist(rescaled_xda)
clf = SpectralClustering(
    n_clusters=4,
    random_state=0,
    gamma=None,
    kmeans_params={"init_max_iter": 5},
    persist_embedding=True,
)

%time clf.fit(X)

/home/runner/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/distributed/client.py:3169: UserWarning: Sending large graph of size 73.62 MiB.
This may cause some slowdown.
Consider scattering data ahead of time and using futures.
  warnings.warn(

/home/runner/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/dask/base.py:1462: UserWarning: Running on a single-machine scheduler when a distributed client is active might lead to unexpected results.
  warnings.warn(

CPU times: user 19.5 s, sys: 2.7 s, total: 22.2 s
Wall time: 27.5 s

SpectralClustering(gamma=None, kmeans_params={'init_max_iter': 5}, n_clusters=4,
                   persist_embedding=True, random_state=0)
In a Jupyter environment, please rerun this cell to show the HTML representation or trust the notebook.
On GitHub, the HTML representation is unable to render, please try loading this page with nbviewer.org.
labels = clf.assign_labels_.labels_.compute()
labels.shape

(1071600,)

labels

array([2, 0, 0, ..., 2, 2, 2], dtype=int32)

Unstack and Visualize

template = flattened_t_xda[:, 0]
output_array = template.copy(data=labels)
unstacked_1988 = output_array.unstack()

/home/runner/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/numpy/core/numeric.py:407: RuntimeWarning: invalid value encountered in cast
  multiarray.copyto(res, fill_value, casting='unsafe')

unstacked_1988

<xarray.DataArray (x: 950, y: 1128)>
array([[2, 0, 0, ..., 0, 0, 0],
       [2, 2, 2, ..., 0, 0, 0],
       [2, 2, 2, ..., 0, 0, 0],
       ...,
       [0, 0, 0, ..., 2, 2, 2],
       [0, 0, 0, ..., 2, 2, 2],
       [0, 0, 0, ..., 2, 2, 2]], dtype=int32)
Coordinates:
  * x            (x) float64 3.353e+05 3.353e+05 ... 3.637e+05 3.638e+05
  * y            (y) float64 4.301e+06 4.301e+06 ... 4.267e+06 4.267e+06
    spatial_ref  int32 32611
    time         datetime64[ns] 1988-07-02T18:03:59.010013
    band         <U3 'red'
Attributes:
    crs:      epsg:32611
raw_plot_1988 = da_1988.sel(band="red").hvplot.image(
    x="x", y="y", geo=True, xlabel="lon", ylabel="lat", datashade=True, cmap="greys", title="Raw 1988"
)

result_plot_1988 = unstacked_1988.hvplot(
    x="x", y="y", cmap="Set3", geo=True, xlabel="lon", ylabel="lat", colorbar=False, title="Spectral Clustering 1988",
)

raw_plot_1988 + result_plot_1988

Spectral Clustering Over Time

Our hypothesis is that the lake’s area is receding over time and so we want to visualize the potential change. Let’s first visually compare the plot of the clustering results from the different time points.

result_plot_1988 + result_plot_2017

By hovering over the lake in the center of each image, we can see that the water was labeled with a cluster label value in both images. Let’s programmatically grab the cluster label at the center for each image.

def get_center_value(arr):
    center_y = arr.shape[0] // 2
    center_x = arr.shape[1] // 2
    center_value = arr[center_y, center_x]
    return int(center_value.values)

water_cluster_1988_label = get_center_value(result_plot_1988.data.value)
water_cluster_2017_label = get_center_value(result_plot_2017.data.value)

Now, we want to align the cluster label for water. Let’s try to set anything that is water to 1 and otherwise 0.

water_1988 = xr.where(unstacked_1988 == water_cluster_1988_label, 1, 0)
water_2017 = xr.where(unstacked_2017 == water_cluster_2017_label, 1, 0)

water_1988.hvplot(geo=True)

WARNING:param.Image02750: Image dimension(s) y and x are not evenly sampled to relative tolerance of 0.001. Please use the QuadMesh element for irregularly sampled data or set a higher tolerance on hv.config.image_rtol or the rtol parameter in the Image constructor.

---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/IPython/core/formatters.py:977, in MimeBundleFormatter.__call__(self, obj, include, exclude)
    974     method = get_real_method(obj, self.print_method)
    976     if method is not None:
--> 977         return method(include=include, exclude=exclude)
    978     return None
    979 else:

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/holoviews/core/dimension.py:1286, in Dimensioned._repr_mimebundle_(self, include, exclude)
   1279 def _repr_mimebundle_(self, include=None, exclude=None):
   1280     """
   1281     Resolves the class hierarchy for the class rendering the
   1282     object using any display hooks registered on Store.display
   1283     hooks.  The output of all registered display_hooks is then
   1284     combined and returned.
   1285     """
-> 1286     return Store.render(self)

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/holoviews/core/options.py:1428, in Store.render(cls, obj)
   1426 data, metadata = {}, {}
   1427 for hook in hooks:
-> 1428     ret = hook(obj)
   1429     if ret is None:
   1430         continue

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/holoviews/ipython/display_hooks.py:287, in pprint_display(obj)
    285 if not ip.display_formatter.formatters['text/plain'].pprint:
    286     return None
--> 287 return display(obj, raw_output=True)

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/holoviews/ipython/display_hooks.py:255, in display(obj, raw_output, **kwargs)
    253 elif isinstance(obj, (CompositeOverlay, ViewableElement)):
    254     with option_state(obj):
--> 255         output = element_display(obj)
    256 elif isinstance(obj, (Layout, NdLayout, AdjointLayout)):
    257     with option_state(obj):

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/holoviews/ipython/display_hooks.py:149, in display_hook.<locals>.wrapped(element)
    147 try:
    148     max_frames = OutputSettings.options['max_frames']
--> 149     mimebundle = fn(element, max_frames=max_frames)
    150     if mimebundle is None:
    151         return {}, {}

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/holoviews/ipython/display_hooks.py:195, in element_display(element, max_frames)
    192 if type(element) not in Store.registry[backend]:
    193     return None
--> 195 return render(element)

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/holoviews/ipython/display_hooks.py:76, in render(obj, **kwargs)
     73 if renderer.fig == 'pdf':
     74     renderer = renderer.instance(fig='png')
---> 76 return renderer.components(obj, **kwargs)

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/holoviews/plotting/renderer.py:396, in Renderer.components(self, obj, fmt, comm, **kwargs)
    393 embed = (not (dynamic or streams or self.widget_mode == 'live') or config.embed)
    395 if embed or config.comms == 'default':
--> 396     return self._render_panel(plot, embed, comm)
    397 return self._render_ipywidget(plot)

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/holoviews/plotting/renderer.py:403, in Renderer._render_panel(self, plot, embed, comm)
    401 doc = Document()
    402 with config.set(embed=embed):
--> 403     model = plot.layout._render_model(doc, comm)
    404 if embed:
    405     return render_model(model, comm)

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/panel/viewable.py:748, in Viewable._render_model(self, doc, comm)
    746 if comm is None:
    747     comm = state._comm_manager.get_server_comm()
--> 748 model = self.get_root(doc, comm)
    750 if self._design and self._design.theme.bokeh_theme:
    751     doc.theme = self._design.theme.bokeh_theme

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/panel/layout/base.py:311, in Panel.get_root(self, doc, comm, preprocess)
    307 def get_root(
    308     self, doc: Optional[Document] = None, comm: Optional[Comm] = None,
    309     preprocess: bool = True
    310 ) -> Model:
--> 311     root = super().get_root(doc, comm, preprocess)
    312     # ALERT: Find a better way to handle this
    313     if hasattr(root, 'styles') and 'overflow-x' in root.styles:

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/panel/viewable.py:666, in Renderable.get_root(self, doc, comm, preprocess)
    664 wrapper = self._design._wrapper(self)
    665 if wrapper is self:
--> 666     root = self._get_model(doc, comm=comm)
    667     if preprocess:
    668         self._preprocess(root)

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/panel/layout/base.py:177, in Panel._get_model(self, doc, root, parent, comm)
    175 root = root or model
    176 self._models[root.ref['id']] = (model, parent)
--> 177 objects, _ = self._get_objects(model, [], doc, root, comm)
    178 props = self._get_properties(doc)
    179 props[self._property_mapping['objects']] = objects

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/panel/layout/base.py:159, in Panel._get_objects(self, model, old_objects, doc, root, comm)
    157 else:
    158     try:
--> 159         child = pane._get_model(doc, root, model, comm)
    160     except RerenderError as e:
    161         if e.layout is not None and e.layout is not self:

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/panel/pane/holoviews.py:423, in HoloViews._get_model(self, doc, root, parent, comm)
    421     plot = self.object
    422 else:
--> 423     plot = self._render(doc, comm, root)
    425 plot.pane = self
    426 backend = plot.renderer.backend

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/panel/pane/holoviews.py:518, in HoloViews._render(self, doc, comm, root)
    515     if comm:
    516         kwargs['comm'] = comm
--> 518 return renderer.get_plot(self.object, **kwargs)

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/holoviews/plotting/bokeh/renderer.py:68, in BokehRenderer.get_plot(self_or_cls, obj, doc, renderer, **kwargs)
     61 @bothmethod
     62 def get_plot(self_or_cls, obj, doc=None, renderer=None, **kwargs):
     63     """
     64     Given a HoloViews Viewable return a corresponding plot instance.
     65     Allows supplying a document attach the plot to, useful when
     66     combining the bokeh model with another plot.
     67     """
---> 68     plot = super().get_plot(obj, doc, renderer, **kwargs)
     69     if plot.document is None:
     70         plot.document = Document() if self_or_cls.notebook_context else curdoc()

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/holoviews/plotting/renderer.py:240, in Renderer.get_plot(self_or_cls, obj, doc, renderer, comm, **kwargs)
    237     defaults = [kd.default for kd in plot.dimensions]
    238     init_key = tuple(v if d is None else d for v, d in
    239                      zip(plot.keys[0], defaults))
--> 240     plot.update(init_key)
    241 else:
    242     plot = obj

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/holoviews/plotting/plot.py:955, in DimensionedPlot.update(self, key)
    953 def update(self, key):
    954     if len(self) == 1 and key in (0, self.keys[0]) and not self.drawn:
--> 955         return self.initialize_plot()
    956     item = self.__getitem__(key)
    957     self.traverse(lambda x: setattr(x, '_updated', True))

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/geoviews/plotting/bokeh/plot.py:108, in GeoPlot.initialize_plot(self, ranges, plot, plots, source)
    106 def initialize_plot(self, ranges=None, plot=None, plots=None, source=None):
    107     opts = {} if isinstance(self, HvOverlayPlot) else {'source': source}
--> 108     fig = super().initialize_plot(ranges, plot, plots, **opts)
    109     if self.geographic and self.show_bounds and not self.overlaid:
    110         from . import GeoShapePlot

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/holoviews/plotting/bokeh/element.py:1891, in ElementPlot.initialize_plot(self, ranges, plot, plots, source)
   1889 if not self.overlaid:
   1890     self._update_plot(key, plot, style_element)
-> 1891     self._update_ranges(style_element, ranges)
   1893 for cb in self.callbacks:
   1894     cb.initialize()

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/geoviews/plotting/bokeh/plot.py:87, in GeoPlot._update_ranges(self, element, ranges)
     86 def _update_ranges(self, element, ranges):
---> 87     super()._update_ranges(element, ranges)
     88     if not self.geographic:
     89         return

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/holoviews/plotting/bokeh/element.py:1048, in ElementPlot._update_ranges(self, element, ranges)
   1045 y_range = self.handles['y_range']
   1046 plot = self.handles['plot']
-> 1048 self._update_main_ranges(element, x_range, y_range, ranges)
   1050 if self._subcoord_overlaid:
   1051     return

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/holoviews/plotting/bokeh/element.py:1190, in ElementPlot._update_main_ranges(self, element, x_range, y_range, ranges)
   1188 if fixed_height:
   1189     plot.frame_height = height
-> 1190     plot.frame_width = int(height/aspect)
   1191     plot.width, plot.height = None, None
   1192 elif fixed_width:

ValueError: cannot convert float NaN to integer

:Image   [y,x]   (value)

water_1988_plot = water_1988.hvplot(
    x="x", y="y", cmap="greys", geo=True, colorbar=False, title="1988 Water"
)

water_2017_plot = water_2017.hvplot(
    x="x", y="y", cmap="greys", geo=True, colorbar=False, title="2017 Water"
)

water_1988_plot + water_2017_plot

WARNING:param.Image02842: Image dimension(s) x and y are not evenly sampled to relative tolerance of 0.001. Please use the QuadMesh element for irregularly sampled data or set a higher tolerance on hv.config.image_rtol or the rtol parameter in the Image constructor.

---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/IPython/core/formatters.py:977, in MimeBundleFormatter.__call__(self, obj, include, exclude)
    974     method = get_real_method(obj, self.print_method)
    976     if method is not None:
--> 977         return method(include=include, exclude=exclude)
    978     return None
    979 else:

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/holoviews/core/dimension.py:1286, in Dimensioned._repr_mimebundle_(self, include, exclude)
   1279 def _repr_mimebundle_(self, include=None, exclude=None):
   1280     """
   1281     Resolves the class hierarchy for the class rendering the
   1282     object using any display hooks registered on Store.display
   1283     hooks.  The output of all registered display_hooks is then
   1284     combined and returned.
   1285     """
-> 1286     return Store.render(self)

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/holoviews/core/options.py:1428, in Store.render(cls, obj)
   1426 data, metadata = {}, {}
   1427 for hook in hooks:
-> 1428     ret = hook(obj)
   1429     if ret is None:
   1430         continue

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/holoviews/ipython/display_hooks.py:287, in pprint_display(obj)
    285 if not ip.display_formatter.formatters['text/plain'].pprint:
    286     return None
--> 287 return display(obj, raw_output=True)

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/holoviews/ipython/display_hooks.py:258, in display(obj, raw_output, **kwargs)
    256 elif isinstance(obj, (Layout, NdLayout, AdjointLayout)):
    257     with option_state(obj):
--> 258         output = layout_display(obj)
    259 elif isinstance(obj, (HoloMap, DynamicMap)):
    260     with option_state(obj):

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/holoviews/ipython/display_hooks.py:149, in display_hook.<locals>.wrapped(element)
    147 try:
    148     max_frames = OutputSettings.options['max_frames']
--> 149     mimebundle = fn(element, max_frames=max_frames)
    150     if mimebundle is None:
    151         return {}, {}

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/holoviews/ipython/display_hooks.py:223, in layout_display(layout, max_frames)
    220     max_frame_warning(max_frames)
    221     return None
--> 223 return render(layout)

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/holoviews/ipython/display_hooks.py:76, in render(obj, **kwargs)
     73 if renderer.fig == 'pdf':
     74     renderer = renderer.instance(fig='png')
---> 76 return renderer.components(obj, **kwargs)

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/holoviews/plotting/renderer.py:396, in Renderer.components(self, obj, fmt, comm, **kwargs)
    393 embed = (not (dynamic or streams or self.widget_mode == 'live') or config.embed)
    395 if embed or config.comms == 'default':
--> 396     return self._render_panel(plot, embed, comm)
    397 return self._render_ipywidget(plot)

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/holoviews/plotting/renderer.py:403, in Renderer._render_panel(self, plot, embed, comm)
    401 doc = Document()
    402 with config.set(embed=embed):
--> 403     model = plot.layout._render_model(doc, comm)
    404 if embed:
    405     return render_model(model, comm)

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/panel/viewable.py:748, in Viewable._render_model(self, doc, comm)
    746 if comm is None:
    747     comm = state._comm_manager.get_server_comm()
--> 748 model = self.get_root(doc, comm)
    750 if self._design and self._design.theme.bokeh_theme:
    751     doc.theme = self._design.theme.bokeh_theme

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/panel/layout/base.py:311, in Panel.get_root(self, doc, comm, preprocess)
    307 def get_root(
    308     self, doc: Optional[Document] = None, comm: Optional[Comm] = None,
    309     preprocess: bool = True
    310 ) -> Model:
--> 311     root = super().get_root(doc, comm, preprocess)
    312     # ALERT: Find a better way to handle this
    313     if hasattr(root, 'styles') and 'overflow-x' in root.styles:

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/panel/viewable.py:666, in Renderable.get_root(self, doc, comm, preprocess)
    664 wrapper = self._design._wrapper(self)
    665 if wrapper is self:
--> 666     root = self._get_model(doc, comm=comm)
    667     if preprocess:
    668         self._preprocess(root)

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/panel/layout/base.py:177, in Panel._get_model(self, doc, root, parent, comm)
    175 root = root or model
    176 self._models[root.ref['id']] = (model, parent)
--> 177 objects, _ = self._get_objects(model, [], doc, root, comm)
    178 props = self._get_properties(doc)
    179 props[self._property_mapping['objects']] = objects

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/panel/layout/base.py:159, in Panel._get_objects(self, model, old_objects, doc, root, comm)
    157 else:
    158     try:
--> 159         child = pane._get_model(doc, root, model, comm)
    160     except RerenderError as e:
    161         if e.layout is not None and e.layout is not self:

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/panel/pane/holoviews.py:423, in HoloViews._get_model(self, doc, root, parent, comm)
    421     plot = self.object
    422 else:
--> 423     plot = self._render(doc, comm, root)
    425 plot.pane = self
    426 backend = plot.renderer.backend

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/panel/pane/holoviews.py:518, in HoloViews._render(self, doc, comm, root)
    515     if comm:
    516         kwargs['comm'] = comm
--> 518 return renderer.get_plot(self.object, **kwargs)

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/holoviews/plotting/bokeh/renderer.py:68, in BokehRenderer.get_plot(self_or_cls, obj, doc, renderer, **kwargs)
     61 @bothmethod
     62 def get_plot(self_or_cls, obj, doc=None, renderer=None, **kwargs):
     63     """
     64     Given a HoloViews Viewable return a corresponding plot instance.
     65     Allows supplying a document attach the plot to, useful when
     66     combining the bokeh model with another plot.
     67     """
---> 68     plot = super().get_plot(obj, doc, renderer, **kwargs)
     69     if plot.document is None:
     70         plot.document = Document() if self_or_cls.notebook_context else curdoc()

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/holoviews/plotting/renderer.py:240, in Renderer.get_plot(self_or_cls, obj, doc, renderer, comm, **kwargs)
    237     defaults = [kd.default for kd in plot.dimensions]
    238     init_key = tuple(v if d is None else d for v, d in
    239                      zip(plot.keys[0], defaults))
--> 240     plot.update(init_key)
    241 else:
    242     plot = obj

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/holoviews/plotting/plot.py:955, in DimensionedPlot.update(self, key)
    953 def update(self, key):
    954     if len(self) == 1 and key in (0, self.keys[0]) and not self.drawn:
--> 955         return self.initialize_plot()
    956     item = self.__getitem__(key)
    957     self.traverse(lambda x: setattr(x, '_updated', True))

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/holoviews/plotting/bokeh/plot.py:945, in LayoutPlot.initialize_plot(self, plots, ranges)
    942     continue
    944 shared_plots = list(passed_plots) if self.shared_axes else None
--> 945 subplots = subplot.initialize_plot(ranges=ranges, plots=shared_plots)
    946 nsubplots = len(subplots)
    948 modes = {sp.sizing_mode for sp in subplots
    949          if sp.sizing_mode not in (None, 'auto', 'fixed')}

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/holoviews/plotting/bokeh/plot.py:1094, in AdjointLayoutPlot.initialize_plot(self, ranges, plots)
   1092     else:
   1093         passed_plots = plots + adjoined_plots
-> 1094         adjoined_plots.append(subplot.initialize_plot(ranges=ranges, plots=passed_plots))
   1095 self.drawn = True
   1096 if not adjoined_plots: adjoined_plots = [None]

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/geoviews/plotting/bokeh/plot.py:108, in GeoPlot.initialize_plot(self, ranges, plot, plots, source)
    106 def initialize_plot(self, ranges=None, plot=None, plots=None, source=None):
    107     opts = {} if isinstance(self, HvOverlayPlot) else {'source': source}
--> 108     fig = super().initialize_plot(ranges, plot, plots, **opts)
    109     if self.geographic and self.show_bounds and not self.overlaid:
    110         from . import GeoShapePlot

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/holoviews/plotting/bokeh/element.py:1891, in ElementPlot.initialize_plot(self, ranges, plot, plots, source)
   1889 if not self.overlaid:
   1890     self._update_plot(key, plot, style_element)
-> 1891     self._update_ranges(style_element, ranges)
   1893 for cb in self.callbacks:
   1894     cb.initialize()

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/geoviews/plotting/bokeh/plot.py:87, in GeoPlot._update_ranges(self, element, ranges)
     86 def _update_ranges(self, element, ranges):
---> 87     super()._update_ranges(element, ranges)
     88     if not self.geographic:
     89         return

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/holoviews/plotting/bokeh/element.py:1048, in ElementPlot._update_ranges(self, element, ranges)
   1045 y_range = self.handles['y_range']
   1046 plot = self.handles['plot']
-> 1048 self._update_main_ranges(element, x_range, y_range, ranges)
   1050 if self._subcoord_overlaid:
   1051     return

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/holoviews/plotting/bokeh/element.py:1190, in ElementPlot._update_main_ranges(self, element, x_range, y_range, ranges)
   1188 if fixed_height:
   1189     plot.frame_height = height
-> 1190     plot.frame_width = int(height/aspect)
   1191     plot.width, plot.height = None, None
   1192 elif fixed_width:

ValueError: cannot convert float NaN to integer

:Layout
   .Image.I  :Image   [x,y]   (value)
   .Image.II :Image   [x,y]   (value)

Now we can take the difference of these water label arrays to see exactly where the water levels has changed.

with xr.set_options(keep_attrs=True):
    water_diff = water_1988 - water_2017

Red pixels (array value ‘1’) of our image below are where water was lost from 1988 to 2017.

water_diff.hvplot(
    x="x", y="y", cmap='coolwarm', geo=True, xlabel="long", ylabel="lat", colorbar=False, title="Water Change 1988-2017",
)

WARNING:param.Image03040: Image dimension(s) x and y are not evenly sampled to relative tolerance of 0.001. Please use the QuadMesh element for irregularly sampled data or set a higher tolerance on hv.config.image_rtol or the rtol parameter in the Image constructor.

---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/IPython/core/formatters.py:977, in MimeBundleFormatter.__call__(self, obj, include, exclude)
    974     method = get_real_method(obj, self.print_method)
    976     if method is not None:
--> 977         return method(include=include, exclude=exclude)
    978     return None
    979 else:

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/holoviews/core/dimension.py:1286, in Dimensioned._repr_mimebundle_(self, include, exclude)
   1279 def _repr_mimebundle_(self, include=None, exclude=None):
   1280     """
   1281     Resolves the class hierarchy for the class rendering the
   1282     object using any display hooks registered on Store.display
   1283     hooks.  The output of all registered display_hooks is then
   1284     combined and returned.
   1285     """
-> 1286     return Store.render(self)

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/holoviews/core/options.py:1428, in Store.render(cls, obj)
   1426 data, metadata = {}, {}
   1427 for hook in hooks:
-> 1428     ret = hook(obj)
   1429     if ret is None:
   1430         continue

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/holoviews/ipython/display_hooks.py:287, in pprint_display(obj)
    285 if not ip.display_formatter.formatters['text/plain'].pprint:
    286     return None
--> 287 return display(obj, raw_output=True)

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/holoviews/ipython/display_hooks.py:255, in display(obj, raw_output, **kwargs)
    253 elif isinstance(obj, (CompositeOverlay, ViewableElement)):
    254     with option_state(obj):
--> 255         output = element_display(obj)
    256 elif isinstance(obj, (Layout, NdLayout, AdjointLayout)):
    257     with option_state(obj):

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/holoviews/ipython/display_hooks.py:149, in display_hook.<locals>.wrapped(element)
    147 try:
    148     max_frames = OutputSettings.options['max_frames']
--> 149     mimebundle = fn(element, max_frames=max_frames)
    150     if mimebundle is None:
    151         return {}, {}

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/holoviews/ipython/display_hooks.py:195, in element_display(element, max_frames)
    192 if type(element) not in Store.registry[backend]:
    193     return None
--> 195 return render(element)

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/holoviews/ipython/display_hooks.py:76, in render(obj, **kwargs)
     73 if renderer.fig == 'pdf':
     74     renderer = renderer.instance(fig='png')
---> 76 return renderer.components(obj, **kwargs)

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/holoviews/plotting/renderer.py:396, in Renderer.components(self, obj, fmt, comm, **kwargs)
    393 embed = (not (dynamic or streams or self.widget_mode == 'live') or config.embed)
    395 if embed or config.comms == 'default':
--> 396     return self._render_panel(plot, embed, comm)
    397 return self._render_ipywidget(plot)

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/holoviews/plotting/renderer.py:403, in Renderer._render_panel(self, plot, embed, comm)
    401 doc = Document()
    402 with config.set(embed=embed):
--> 403     model = plot.layout._render_model(doc, comm)
    404 if embed:
    405     return render_model(model, comm)

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/panel/viewable.py:748, in Viewable._render_model(self, doc, comm)
    746 if comm is None:
    747     comm = state._comm_manager.get_server_comm()
--> 748 model = self.get_root(doc, comm)
    750 if self._design and self._design.theme.bokeh_theme:
    751     doc.theme = self._design.theme.bokeh_theme

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/panel/layout/base.py:311, in Panel.get_root(self, doc, comm, preprocess)
    307 def get_root(
    308     self, doc: Optional[Document] = None, comm: Optional[Comm] = None,
    309     preprocess: bool = True
    310 ) -> Model:
--> 311     root = super().get_root(doc, comm, preprocess)
    312     # ALERT: Find a better way to handle this
    313     if hasattr(root, 'styles') and 'overflow-x' in root.styles:

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/panel/viewable.py:666, in Renderable.get_root(self, doc, comm, preprocess)
    664 wrapper = self._design._wrapper(self)
    665 if wrapper is self:
--> 666     root = self._get_model(doc, comm=comm)
    667     if preprocess:
    668         self._preprocess(root)

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/panel/layout/base.py:177, in Panel._get_model(self, doc, root, parent, comm)
    175 root = root or model
    176 self._models[root.ref['id']] = (model, parent)
--> 177 objects, _ = self._get_objects(model, [], doc, root, comm)
    178 props = self._get_properties(doc)
    179 props[self._property_mapping['objects']] = objects

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/panel/layout/base.py:159, in Panel._get_objects(self, model, old_objects, doc, root, comm)
    157 else:
    158     try:
--> 159         child = pane._get_model(doc, root, model, comm)
    160     except RerenderError as e:
    161         if e.layout is not None and e.layout is not self:

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/panel/pane/holoviews.py:423, in HoloViews._get_model(self, doc, root, parent, comm)
    421     plot = self.object
    422 else:
--> 423     plot = self._render(doc, comm, root)
    425 plot.pane = self
    426 backend = plot.renderer.backend

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/panel/pane/holoviews.py:518, in HoloViews._render(self, doc, comm, root)
    515     if comm:
    516         kwargs['comm'] = comm
--> 518 return renderer.get_plot(self.object, **kwargs)

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/holoviews/plotting/bokeh/renderer.py:68, in BokehRenderer.get_plot(self_or_cls, obj, doc, renderer, **kwargs)
     61 @bothmethod
     62 def get_plot(self_or_cls, obj, doc=None, renderer=None, **kwargs):
     63     """
     64     Given a HoloViews Viewable return a corresponding plot instance.
     65     Allows supplying a document attach the plot to, useful when
     66     combining the bokeh model with another plot.
     67     """
---> 68     plot = super().get_plot(obj, doc, renderer, **kwargs)
     69     if plot.document is None:
     70         plot.document = Document() if self_or_cls.notebook_context else curdoc()

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/holoviews/plotting/renderer.py:240, in Renderer.get_plot(self_or_cls, obj, doc, renderer, comm, **kwargs)
    237     defaults = [kd.default for kd in plot.dimensions]
    238     init_key = tuple(v if d is None else d for v, d in
    239                      zip(plot.keys[0], defaults))
--> 240     plot.update(init_key)
    241 else:
    242     plot = obj

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/holoviews/plotting/plot.py:955, in DimensionedPlot.update(self, key)
    953 def update(self, key):
    954     if len(self) == 1 and key in (0, self.keys[0]) and not self.drawn:
--> 955         return self.initialize_plot()
    956     item = self.__getitem__(key)
    957     self.traverse(lambda x: setattr(x, '_updated', True))

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/geoviews/plotting/bokeh/plot.py:108, in GeoPlot.initialize_plot(self, ranges, plot, plots, source)
    106 def initialize_plot(self, ranges=None, plot=None, plots=None, source=None):
    107     opts = {} if isinstance(self, HvOverlayPlot) else {'source': source}
--> 108     fig = super().initialize_plot(ranges, plot, plots, **opts)
    109     if self.geographic and self.show_bounds and not self.overlaid:
    110         from . import GeoShapePlot

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/holoviews/plotting/bokeh/element.py:1891, in ElementPlot.initialize_plot(self, ranges, plot, plots, source)
   1889 if not self.overlaid:
   1890     self._update_plot(key, plot, style_element)
-> 1891     self._update_ranges(style_element, ranges)
   1893 for cb in self.callbacks:
   1894     cb.initialize()

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/geoviews/plotting/bokeh/plot.py:87, in GeoPlot._update_ranges(self, element, ranges)
     86 def _update_ranges(self, element, ranges):
---> 87     super()._update_ranges(element, ranges)
     88     if not self.geographic:
     89         return

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/holoviews/plotting/bokeh/element.py:1048, in ElementPlot._update_ranges(self, element, ranges)
   1045 y_range = self.handles['y_range']
   1046 plot = self.handles['plot']
-> 1048 self._update_main_ranges(element, x_range, y_range, ranges)
   1050 if self._subcoord_overlaid:
   1051     return

File ~/miniconda3/envs/cookbook-dev/lib/python3.10/site-packages/holoviews/plotting/bokeh/element.py:1190, in ElementPlot._update_main_ranges(self, element, x_range, y_range, ranges)
   1188 if fixed_height:
   1189     plot.frame_height = height
-> 1190     plot.frame_width = int(height/aspect)
   1191     plot.width, plot.height = None, None
   1192 elif fixed_width:

ValueError: cannot convert float NaN to integer

:Image   [x,y]   (value)

We did it! We are observing the change in the lake shoreline over time using a simple spectral clustering approach.

Let’s finish things off by adding some geo tiles as a background. To only display the colored pixels overlaid on geo tiles, we could either set the array’s background value (‘0’) to ‘Not a Number’ (NaN), or we could just inform hvPlot that we want the background valued pixels to be transparent with .redim.nodata(value=0).

water_diff.hvplot(
        x="x", y="y", width=400, height=400, cmap='coolwarm', geo=True, xlabel="lon", ylabel="lat", alpha=1, colorbar=False, title="Water Loss from 1988 to 2017", tiles="ESRI",
).redim.nodata(value=0)

WARNING:param.Image03132: Image dimension(s) x and y are not evenly sampled to relative tolerance of 0.001. Please use the QuadMesh element for irregularly sampled data or set a higher tolerance on hv.config.image_rtol or the rtol parameter in the Image constructor.

Summary

Starting from raw Landsat data, we have used a simple spectral clustering approach to observe the change in a lake water’s extent across time.

What’s next?

Adapt this notebook for your own use case or select another workflow example notebook.

Resources and References

  • Authored by Demetris Roumis circa Jan, 2023

previous

Preprocessing - Regrid